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is shown that for certain source strengths one class of regular solutions approaches 
Minkowski spacetime, while the other comes arbitrarily close to black hole solutions. 
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1. Introduction 

Static bodies resulting from combined gravitational attraction and electrical 
(electrostatic) repulsion have been playing a central role in many investigations, 
eg. P3 El IHj and references therein, connected with general relativistic treatment of 
gravitational contraction. As a part of Einstein equation, one has to specify energy- 
momentum tensor which, for a perfect fluid, has the form 

T^v = {p + p)u^u v +pg^ u , (1) 

where p, p, and g^ v are mass density, pressure, four-velocity and metric (in 
geometrized units: G — 1 = c). In addition, one has to specify the equation of state 
relating the pressure to the density p = p(p). According to the Oppenheimer-Snyder 
scenario jlj, in a static spherically symmetric perfect fluid ball with the above energy 
momentum tensor, after stationary phase and at the end of nuclear burning, the central 
pressure is no longer able to counterbalance the gravitational attraction and the ball 
starts to collapse. By setting p = the energy momentum tensor assumes the form 
characteristic for dust ball = pu^u v . Eventually, as the ball shrinks one can expect 
formation of a black hole. 

The Majumdar-Papapetrou formalism El 13 IE] can be applied to study 
properties of spacetimes generated by T^ v = pu^u u where p represents static electrically 
counterpoised dust (or extremal charged dust, ECD), i.e. pressureless matter in 
equilibrium under its own gravitational attraction and electrical repulsion 121 El- 
Special attention has been paid to ECD distributions leading to spacetimes that are 
regular everywhere, but with exterior arbitrarily close to that of an extremal Reissner- 
Nordstr0m (ERN) black hole QUI EI]. 

This paper considers the coupled Einstein-Maxwell field equations in the 
Majumdar-Papapetrou approach. With the energy-momentum tensor given by (0) with 
p = 0, our study will concentrate on diverse (analytic) spherically symmetric forms of 
ECD distributions and corresponding solutions will be found. The theoretical framework 
is presented in section |2l In section El we present the static solutions to Einstein- 
Maxwell equations for diverse ECD distributions obtained numerically, and in section 0] 
for the 5-shell distribution, all of which show interesting bifurcating properties. In 
sectional we discuss the bifurcating behaviour and show that diverse distributions can be 
treated on equal footing. Similar behaviour found in gauge theories and Einstein- Yang- 
Mills systems is discussed, and we give some details on numerical methods for finding 
bifurcating solutions. In section EJ we show that with adjusting of the source strength 
and appropriate rescaling of the solutions, the spacetimes with exteriors arbitrarily close 
to the ERN case can be obtained. Conclusions are given in section and extension of 
the present work is proposed. 
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2. Theoretical framework 

One usually starts by denning a static spherically symmetric line element of the form 

ds 2 = -B{r)dt 2 + A{r)dr 2 + r 2 du 2 , (2) 
where da; 2 = d$ 2 + sin 2 $dy? 2 . For an asymptotically flat spacetime, one requires 

^)Uoo = ^(r)U 00 = l, (3) 

and for a nonsingular spacetime A(r = 0) = 1. Papapetrou has shown [UJ [7j that the 
line element can be written as 

ds 2 = -e* dt 2 + e^ (dX 2 + dY 2 + dZ 2 ) (4) 

and that it is possible to connect functions x an d f in such a way as to satisfy (as one 
possibility) dx/d<p = —1 or x(f) — ~ L P- Together with Majumdar's assumption |3] 
about connection between g tt (or g 00 ) metric component and a scalar component of the 
electromagnetic potential comprising the electromagnetic energy-momentum tensor 
T^ m ' ) , the line element ds 2 in the Majumdar-Papapetrou form can be written in the 
harmonic coordinates (t, X) as 

ds 2 = -e" v dt 2 + e v (dX 2 + dY 2 + dZ 2 ) (5) 

(in the original notation of [ZD- We will assume the spherical symmetry as a general 
symmetry requirement for our problem. The essential ingredient to the line element 
(0), as shown by Majumdar, are additional assumptions which should be made about a 
specific form for the energy-momentum tensor that enters the Einstein equations 
1 

R^u - -g^R = SnT^ (6) 

(in geometrized units). The component of T^ u due to electromagnetic fields is given by 

i 

tH — p °~ p n pP° p (7) 

where F^ v = A u ^ — A^^ is the antisymmetric electromagnetic field tensor, comma 
denoting the ordinary derivative to differ from semicolon denoting the covariant 
derivative. It satisfies the empty space equation F^ v v = 0. 

A more general situation differing from the one described above corresponds to 

T tw = T^+T^, (8) 

where the superscript (m) denotes matter. The absence of leads to Reissner- 
Nordstr0m solutions with B(r) = A(r)~ l = (1 — 2m/ r + q 2 /r 2 ) for A = e/r. It is 
shown in Ref. that for charged sphere of radius a, for r large enough, m could be 
set equal to zero, whereas for the interior solution, i.e. r < a, one is not permitted 
to put m = 0. It follows that electric charge contributes to the gravitational mass of 
the system. In addition, mass as given in B(r) = A{r)~ x above satisfies a positivity 
condition m > 0. Citing Ref. [12. 'a charged sphere must have a positive mass'. 

In the case of the complete T^ u given in (jHJ) one assumes also a perfect fluid of 
classical hydrodynamics for which, in general, energy-momentum tensor is given by . 



Regular and quasi black hole solutions for charged dust A 

With p = 0, the energy momentum tensor ([!} reduces to = pu^u v , where p is now 
an invariant dust density and is the four-velocity which will be assumed to be in 
non-moving (co-moving) coordinate system so that 

u» = ^-^u» = 5{;(-g 00 )-^. (9) 
as 

The electromagnetic part of is given, as before, with (jJJ). The four-current density 
of the the charge/matter /dust distribution is given by 

r = Pe u». (io) 

At this point, the crucial assumption is made by setting the electric charge density p e 
equal to the matter density p, or somewhat more general p e = ±p. Now the metric (|5]l 
can be written in the form 

ds 2 = -U- 2 dt 2 + U 2 (dX 2 + dY 2 + dZ 2 ) . (11) 

The Einstein field equations incorporating the above introduced ingredients are given 
by 



in 



\F/F va - //,„-/••""/•;, I + Pelt/," i 



(12) 



and by using the Majumdar (or Majumdar-Papapetrou) assertion, the generalized 
(nonlinear) Poisson equation (with p = p e ) 

F^. v = ATtJ^ = At: pu^ (13) 

assumes a simple form 

V 2 U = -An pU 3 . (14) 

Here U and p are functions of the co-moving coordinates (X, Y, Z) determined from the 
zeroth component of the electromagnetic potential A^, i.e. 

A" = 5Z<t> and U = — *— . (15) 

1 - (f) 

Confining our attention to the spherically symmetric case, the Majumdar-Papapetrou 
metric can be written as 

ds 2 = -U' 2 {R) dt 2 + U 2 (R) (dR 2 + R 2 dVl 2 ) . (16) 

The metric function U is a function of R only, so (|14j) reduces to 

V 2 U(R) = R- 2 (R 2 U'(R))' = -Anp(R)U 3 (R), (17) 

where the prime denotes the differentiation with respect to R. The line element ()16|) 
can be expressed in the standard form (J2J through the coordinate transformation 

r = RU(R), (18) 

with the following relations among the profile functions: 

B(r) = U- 2 (R) (19) 
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and 



1 R dU(R) , , 

1 + -777-- ~r ■ (20) 



y/A{r) U(R) dR 

In regions of space where p(R) = 0, the nonlinear equation (f*T7j) reduces to a 
homogeneous equation with the general solution 

Tfi 

U(R) = k + -, (21) 

where k and m are integration constants. Using (|18j) - (j2(Jj) . to express the line element 
in the standard form (J2J), one obtains B{r) = A(r)~ l = k~ 2 (1 — m/r) 2 . If the region of 
space we are considering extends to r — ► oo, according to the requirements (jHJ), we set 
k — 1, and only m remains as a free parameter. The line element is then 

m\ 2 2 / m\~ 2 



ds 2 = -[}--) dt 2 +(}--) dr + r d ^ > ( 22 ) 

which we recognize as the extremal Reissner-Nordstr0m (ERN) spacetime. The general 
Reissner-Nordstr0m (RN) solution specifies the asymptotically flat spacetime metric 
around an electrically charged source. Expressing the RN line element in the standard 
form (j2J), one has B(r) = A(r)~ l = (1 — 2m/r + q 2 /r 2 ), where m and q are the ADM 
mass and charge of the source. For m < q, the RN spacetime is regular everywhere 
except at r = 0, while for m > q > 0, in addition to the singularity at r = 0, it 
exhibits two horizons located at r = m ± (m 2 — q 2 ) 1 ^ 2 . In the ERN case, i.e. m = q 
(gravitational attraction balances electrical repulsion), only one horizon at r = m is 
present. It is important to note that the solution (|21j) is valid in the range < R < oo 
which is according to eq. (fT£J) mapped to m < r < oo. That is to say that the solution 
(|21[) specifies only the external part of the ERN spacetime (see [TU]). 

We shall proceed to solve the field equation (fTTjl for several distributions of 
electrically counterpoised dust (or extremal charged dust, ECD) p(R). The ECD 
distributions that will be considered effectively vanish at large R, so we expect our 
spacetimes to behave like for R — >• oo. The parameter moo = R(U — 1) |_r^oo, 
characterizing the asymptotic behaviour of U(R) at large R, is the ADM mass seen by 
the distant observer. The behaviour of U(R) as R — > can be characterized by the 
parameter m = RU \r-^o- The metric (|T5|) is regular at R = only if mo = and, 
according to (jl8j) . the range < R < oo is mapped to < r < oo. On the other 
hand, if mo > the function U(R) is infinite at R = which indicates the singularity. 
As the range < R < oo is now mapped to m < r < oo, this singularity is, in the 
r-coordinate, located at r = mo- Therefore, we may understand the parameter mo as 
the 'mass below horizon'. If p = we have m = m^, while in case of non- vanishing 
ECD density p we have 

moo = m + m p , (23) 

where m p can be understood as the contribution of ECD to the ADM mass of the 
configuration. The 'ECD mass' m p is the space integral of the ECD density p: 

poo roo 

4vr / pU 3 R 2 dR = 4n pA 1/2 r 2 dr. (24) 

J J ran 
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We conclude this Section by pointing out, for later convenience, that if certain U(R) 
and p(R) solve the equation (fTTj) . one is allowed to rescale the functions 

U(R) — > U*(R) = U(aR), 

p(R) — > p*(R) = a 2 p(aR), 1 } 

which, as a consequence, rescale the mass parameters according to 

m — ► in* = m/a. (26) 

These relations will help us compare field configurations corresponding to different mass 
parameters. 



3. General ECD distributions and bifurcating solutions 

Our approach in constructing solutions to the Majumdar-Papapetrou systems is to 
assume certain ECD distribution p(R) and then integrate the differential equation (fTTj) 
numerically to obtain the metric function U(R), and is therefore different from the 
approach used in [TJ [TQj where p(R) is analytically reconstructed from the assumed 
metric and the requirement that p(R) may not be negative. At the expense of requiring 
numerical procedures our approach allows complete freedom in choosing the shape of 
p(R), which we found more natural. However, in accord with the requirement that 
our solutions be asymptotically flat, we choose ECD distributions p(R) that are well 
localized in the -R-coordinate. More specifically, we require that the total mass/charge 
in a linear theory 47r J °° pR 2 dR is finite, i.e. that p(R) falls off more rapidly than 1/R 3 . 
As our first example we take p{R) of the form 

P(R) = (R/R) e _JVA , (27) 

where R is the length scale that in further text we set equal to 1, rj has the role of 
an adjustable source strength factor (in appropriate units) and 2Att is a normalization 
constant that renders An f °° pR 2 dR = 1 for n = 1 and R — 1. 

As (11 7J) is a second order differential equation we have to impose two boundary 
conditions onto a solution. As the first boundary condition, we use the requirement 
that the spacetime is asymptotically flat, i.e. U \r->oo — 1> an d as the second we chose 
to fix the parameter mo = RU |_r^o- We obtain the solutions by numerical integration 
[T3] of the differential equation. Starting with the solution of (j!7)) for source strength 
rj = 0, i.e. moo = m o, we slowly increase the value of rj. Increasing of rj evolves the 
solution toward m^ > mo, but one may only increase the value of t] up to a critical 
value r] c . The dependence of m^ on i] in solutions obtained in this way is shown as the 
lower part of the solution tracks in figure H The critical points are labeled (a), (d), (g), 
and (h) , and the corresponding values of rj c and moo are given in table ^ 

In addition to these solutions there exists another class of solutions: starting from 
a critical solution, i.e. one obtained with r] c , in our numerical procedure, we could 
require a slight increase in moo and allow the source strength r] to adjust itself freely. 
The solutions to the differential equation can be found and it turns out that, for this 
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Figure 1. The dependence of moo (ADM mass) on source strength 77 in solutions to 
(fT7|) with ECD distribution given by l|T7J) and boundary conditions m = (leading to 
regular spacetimes) and mo = 1,2,4 (leading to singular spacetimes). The numerical 
values of 77 and for the solutions labeled by letters (a)-(h) are given in tabled 



Table 1. Numerical values of 77 and m^ for the solutions labeled by letters in figure^ 



m 


Vc 


TOoo(?? c ) 









0.589 272 


1.80105 (a) 


0.795196 (6) 


3.713 25 (c) 


1 


0.265 212 


3.266 61 (d) 


1.976 69 (e) 


5.849 51 (/) 


2 


0.137 966 


4.906 04 (g) 






4 


0.053 701 


8.394 04 (h) 







class of solutions, as we increase moo, the source strength 77 decreases! These solutions 
are indicated as the upper part of the tracks in the versus 77 diagram (figure Q). 
Therefore, for source strength 77 > r] c there appears to be no solution to (fT7|) that is 
asymptotically flat, while if starting from rj c (i.e., from the critical solution such as (a) 
or (d) in figure Q]), the solutions bifurcate by following either lower (through points (b) 
or (e)) or upper moo-branch (through points (c) or (/)). This kind of bifurcation is of 
the 'turning point type' as explained in jHj (see also [To]). 

In figure ^ the critical point on the curve corresponding to the boundary condition 
m = is obtained with rj c = 0.589 27 and is labeled (a). The points (6) on the lower 
and (c) on the upper moo-branch indicate the two independent solutions obtained with 
77 = |7/ c = 0.47141. The corresponding masses are given in tabled The components 
of the metric in the r-coordinate A(r) = g rr and B(r) = g u for the solutions labeled 
(a)-(c), are shown in figure|21 At large r, the metric components A(r) and B(r) coalesce 
into the ERN metric (|22j) . But in contrast to the ERN metric, A(r) and B(r) of our 
solutions are finite at all r so the spacetime does not involve an event horizon. This 
is true for all solutions obtained with the boundary condition mo = 0. We call them 
the regular solutions, although in section |H1 we show that these solutions may come 
arbitrarily close to having an extremal horizon. 
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r 



1 5 10 



Figure 2. Regular spacetime metric components g~} — A~ Y (r) (thick lines) and 
gtt — B{r) (dashed lines) obtained by solving the field equation (|17fl with the ECD 
distribution given by lt?f|) and boundary condition m = 0. The solution (a) is the 
critical solution, while the solutions (b) and (c) are the two independent solutions 
obtained with r\ = |?7 C (see figure ^ an( l table for numerical values). The ERN 
metric components A^ 1 = B = (1 — rtioo/r) 2 are shown for corresponding moo values 
(dotted lines). 



B, A~ l 




r 



1 5 10 



Figure 3. Regular spacetime metric components — A~ l (r) (thick lines) and 
g tt — B(r) (dashed lines) obtained by solving the field equation (|17|l with the ECD 
distribution given by l|27|) and boundary condition mo = 1. The solution (d) is the 
critical solution, while the solutions (e) and (/) are the two independent solutions are 
obtained for 77 = ^rj c , (see figure^ an d table^for numerical values). The ERN profile 
functions A^ 1 = B = (1 — m^/r) 2 are shown for corresponding moo values (dotted 
lines). 
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The situation is substantially different for the solutions obtained with the boundary 
condition m > 0. The metric components for the case mo = 1, labeled (d)-(f), in 
figure ^ and tabled are shown in figure El At r = mo, the metric component A(r) = g rr 
diverges (1/A(r) reaches zero) which indicates an event horizon. When approached from 
r > mo the functions 1/A(r) and B(r) appear to have a double zero at r = and in this 
sense this spacetime singularity is equivalent to the ERN event horizon. Recall that at 
r < m the spacetime metric is not specified by the solution to the field equation (JT7J). 
This type of solutions we call singular. 

For different choices of ECD distribution p, as for instance p(R) oc R 2 exp(— R 2 ), 
basically the same bifurcating behaviour of solutions occurs. This fact will be used to 
simplify the treatment of p's in section 03 



4. The 5-shell ECD distribution 

Here we will consider the case of ECD distributed on a spherical shell of radius _R - A 
situation similar to this one was considered in jTHI, while thick shells were considered in 
HI]. We set 

p(R)=r ] 5(R-R ). (28) 

Both in the interior (R < Rq) and in the exterior (R > Rq) space the general solution 
to (JT7J) is of the form (|2~Tj). We set the exterior solution to be asymptotically flat and 
allow the interior solution to be singular at R = 0: 

TTim _ J k + m /R = U h R<R , 

u w-\ 1 + moo/R = u E , R>Ro . W 

The requirement that U\ — C/e at R — Rq fixes the value of k: 

k = 1 + (rrioo - m )/R . (30) 

Integration of the differential equation (fT7|) with the rhs involving the 5-shell source (|2*H|) 
yields 

R 2 U ^R+e^_ A7tR 2 r]U ^3 (31) 

When the solution (J2T?)) is substituted into the above relation the source strength rj, m , 
and the position of the 5-shell source, are interrelated: 

moo ~ m ( moo V 3 /QO s 

This relation leads to bifurcating solutions similar to those discussed in section [3 

We now restrict the discussion to the regular solutions, i.e., we set m = and 
proceed with only one mass parameter m = m p = given by equation (J24j) . The 
source strength is 



"=4^(TT^ (33) 
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1 



21 n R 




in 



1/2 



Figure 4. The source strength rj vs. £ 
equation l(T7|) with 5-shell source <|28|) . 



m/Ro in regular solutions to the field 



where £ = m/R and it follows that at £ c = 1/2 there is a maximum r^ c = (277r.Ro) 1 as 
shown in figure 0] For the regular solution the .R-coordinate metric component U(R) is 



According to (|18jl. the position of the 5-shell in the r-space is ro = m + Rq. In the 
interior of the 5-shell the r-space metric components are constants 



while in the exterior the metric follows the ERN metric given by (|22|) . The metric 
components of only the critical regular solution involving the 5-shell ECD distribution 
are shown in figure El 

5. Bifurcating behaviour and unified treatment of ECD distributions 

Solutions for different ECD distributions p(R) found in sections El and |U exhibit 
bifurcating behaviour with respect to the source strength parameter rj. When rj < rj c , one 
is able to find two independent /different solutions to the same (nonlinear) differential 
equation with the same boundary conditions. As rj approaches the critical value r] c the 
two solutions become identical. While the asymptotic flatness of the two independent 
solutions is fixed by the boundary conditions, their ADM mass ()24|) is different. 

If one takes any versus rj curve in figure Q (and corresponding spacetime metric 
component solutions), it is not possible to distinguish the upper branch from the lower 
branch solutions apart from calculated from (|24j) . Only a sequence of solutions 
provides a means to determine the branch along which solutions propagate. 

One can expect that metric component configurations leading to smaller ADM 
mass when a mass source strength 7] is decreased are physically more acceptable then 




(34) 




(35) 
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1 





r 



m 



r = Ro+m 



Figure 5. The metric components g~ r x = A _1 (r) (thick line) and gu = B{r) (dashed 
line) for the regular critical solution of the field equation l(T7l) with <5-shell source . 
In the exterior of the shell located at ro = m + Ro, the metric is equivalent to the ERN 
(E2> for the mass m (extended by dotted line into the interior space). 



their counterparts which produce larger and larger mass = ttiadm- Indeed, when 
the source strength is decreased from its critical value the low-mass bifurcation branch 
leads to classical (Newtonian gravity) configurations of ECD. One may say that the 
upper branch solutions (points (c), (/), (h) in figure 1) may tend to infinite mass 
which renders them to be (physically) unacceptable. Further investigation of stability 
properties of these solutions might give some definite answer to the above assertions. 
On the other hand, the spacetimes corresponding to high-mass configurations have the 
interesting and important property of being quasi-singular which we discuss further in 
section |UJ 

Based on the bifurcation properties of the solutions to the field equation (fTTj) . we 
are going to introduce a normalization of the ECD distributions appropriate for the 
Majumdar-Papapetrou (MP) formalism. This will allow us to treat the diverse ECD 
distributions on equal footing and more easily explore their properties. We will focus 
only on the regular solutions, so only one mass parameter m = m p = will be used. 
The ECD distributions considered in section El will be used to generate spacetimes that 
come arbitrarily close to having extremal horizons. 

We first consider a general distribution p(R) = i] p (R) where r\ is the source 
strength factor and po is normalized so that 



In the linear theory, the contribution of p = r] p to the total mass and charge of the 
configuration would be m = q = r], and there would be no upper bound imposed onto 
the source strength r\ (see figure dotted line). In our case where the field equation 
(fTTjl is non-linear, solutions exist only for values of i] less or equal to a critical value 
r/ c . For r] = r] c the solution is unique, we call it the critical solution, and we label 




(36) 
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0.2 0.4 0.6 



Figure 6. Dependence of mass m on the source strength r\ for the regular solutions 
to the field equation l|17|) with ECD distribution given by l|37|) . Going from left to 
right, the tracks are ordered as in the table critical solutions are indicated by circles. 
Dotted line indicates the expected behaviour in the linear theory. 



the corresponding mass with the symbol m c . For 7] < 7] c there are two independent, 
bifurcating solutions, leading to masses that we label where it holds m_ < m c < m + . 
As examples of ECD distributions normalized according to (}36j) . we considered 

Po(R) = n(a, b) (R/R) a exp (-(/?/ R) b ), (37) 

where as before R = 1 (i.e. a unit mass), a, b > are parameters and n(a,b) = 
bj (47iT[(3 + a)/b]). The dependence of m on t] for the parameter values a = 0, 1, 2 and 
b = 1,2 is shown in figure El The numerical values for r] c and m c are given in tabled 
As can be seen, only the m~-branch, and only at the low source strength limit, behaves 
as it would in the linear theory. 

In the context of the field equation ()17p . we can formulate a more natural approach 
to normalization of the sources for the Majumdar-Papapetrou spacetimes that we will 
call MP-normalisation. For the MP- normalized ECD distribution p(R), we require that 
the critical source strength factor and corresponding mass are both equal to unity, i.e. 
T] c = m c = 1. To obtain a MP-normalized ECD distribution starting from an arbitrarily 
normalized ECD distribution po(R), we first solve (|T7|) to obtain the values i] c and m c . 
The MP-normalized source p is then constructed by rescaling the density po according 
to (US)) and ()26|) with the scaling parameter a = m c . The MP-normalized source p 
corresponding to po becomes 

p(R) = m 2 c T] c p (m c R) . (38) 

As an example we can take the 5-shell ECD distribution ()28|) discussed in section 0] for 
which we obtained r\ c = (277r_R ) _1 and m c = Ro/2 (see figure EJ). The MP-normalized 
5-shell density follows as p{R) = (l/547r) S(R — 2). In the case of the ECD distributions 
()37|) . the MP-normalized version reads 

p(R) = m a c +2 r] c n{a, b) R a exp(-(m c i?) 6 ), (39) 
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Table 2. The critical points in the regular solutions to (|17|) obtained with the with 
the ECD distribution ijTfjl. 



a 


b 


l/n(a, b) 




m c 





2 


7T 3 / 2 


0.173 340 


0.542 345 


1 


2 


2tt 


0.209 013 


0.670181 


2 


2 


3^3/2 
2 " 


0.238 647 


0.774 735 





1 


8tt 


0.416177 


1.21882 


1 


1 


24tt 


0.589 272 


1.80105 


2 


1 


96tt 


0.758 753 


2.37119 



where r] c and m c obtained numerically are given in table [2] 

The bifurcating solutions found in sections H3 and 0] have similar properties to the 
solutions found long time ago in another nonlinear field theory, i.e. in Yang-Mills gauge 
theory with external sources [X8j. Stability of solutions has been investigated and 
different stability properties with respect to radial oscillations has been found I2U]. 
In the context of Einstein- Yang-Mills theory some of solutions have been found in the 
gravitating SU(2) monopole context [2H 122 I2H 123] and an SU(3) extension of this 
model |2S], which could be characterized as bifurcating. Similar behaviour has recently 
been discussed in an extension of the Standard Model [2H] • In some of the above papers 
the term bifurcation has been used rather loosely because some solutions found there 
do not follow the requirements from literature fSl IHj- Although an explanation of 
such behaviour is missing in all those works, a careful numerical analysis done in [21] as 
well as [27] offers an opportunity to compare solutions discussed before and bifurcating 
behaviour of our solutions. 

We end this section by giving some technical details related to the numerical 
procedure we used to generate the bifurcating solutions: assume a sequence of solutions 
fi(x;r]i), % = l,2..c to a nonlinear differential equation obtained with a sequence of 
source strength factors r\i. In order to find an expected upper or lower branch solution, 
when a solution is already found, one can replace the parameter rj by a function t] (x) 
described by differential equation dr](x) /dx = 0. Since the system is now enlarged by one 
first order differential equation, one additional boundary condition has to be supplied. 
This boundary condition can be constructed from the functional value fk(x s ,r]k), by 
requiring that f(x s ,r](x s )) = fk(x s ,r]k) ± Sf, where f(x,rj(x)) is a part of the enlarged 
system, and Sf (being a small quantity) produces upper (+) or lower (— ) functional 
value, giving an upper or lower branch solution. It is reasonable to choose fk(x s ,r]k) 
close to the critical value r] c . Once a solution is found the simple continuation used in 
the numerical code [13] could produce a corresponding sequence of solutions. Relative 
tolerances could be chosen to be very stringent (typical order is 10 -6 to 10 -8 ) which 
assures that bifurcation is not produced by loose numerical boundaries and justifies the 
six significant digits given in table [21 
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~~ T~ i 2 3 4 5 

Figure 7. The r-space metric components 1/A and -B (dashed lines) for the critical 
regular solutions obtained with the first three MP-normalized sources of figure H3 and 
table [5] (dashed lines), and for m = 1 ERN metric components (dotted line). (Also to 
be compared with metric components in figure 0) 

6. ECD distributions and quasi black holes 

Using a MP-normalized ECD distribution p, the critical regular solution to the field 
equation (jT7j) is obtained with the source strength 77 = t] c = 1. The r-coordinate metric 
components of the critical regular solutions to (JT7J) and the ECD distribution (}39|) with 
parameters a = 2 and b = 0, 1,2 are shown in figure [7| At large r, metric components 
A and B coalesce into the m — 1 ERN metric, while at the intermediate values where 
major part of the ECD is distributed, they are manifestly regular and do not show 
significant dependence on the choice of the shape of the ECD distribution. By using 
77 < 1, these solutions can, due to the bifurcating behaviour, either evolve toward flat 
space along m~ branch, or toward higher mass configurations along the m + branch. Any 
of these solutions can be rescaled to describe a unit mass configuration if the rescaling 
according to (|25jl and (|2fij1 is carried out with the scaling parameter a = m. Starting 
from the critical solutions shown in figure [7| we followed the m + -branch to obtain the 
field configurations corresponding to m + = 10 and m + = 100 which we then rescaled 
to restore the unit mass configurations. The metric components of these solutions are 
shown in figure |H1 

It is plausible that, by following the m + -branch toward higher masses and then 
rescaling the solutions to describe the m = 1 configurations, the minimum of the function 
1/A(r) is deeper and closer to r = 1. Asymptotically, in the region r > 1 the metric 
components A(r) and B(r) follow the ERN metric and as the point r = 1 is approached 
from r > 1 there appears to be a double zero both in 1/A(r) and B(r). In the r < 1 
region we would asymptotically have B = which would not allow for timelike intervals. 
However, since our spacetimes are regular everywhere and such a situation is realized 
only asymptotically, we can consider a radially moving photon for which ds 2 = and 
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Figure 8. The r-space metric components obtained by using the MP-normalized 
sources as in figure[7|to generate m + — 10 (upper panel) and m + = 100 (lower panel) 
regular solutions which were then rescaled to describe unit mass objects, m = 1 ERN 
metric (dotted line). 



duj 2 = 0, so dt/dr = ±^ A/B. The ratio of the metric components a/ A/B for the 
m + = 10/rescaled solutions, related to the time required for the photon to transverse a 
unit distance in the r-coordinate, is shown in figure H3 

The distributions of the ECD for the solutions of figure |H1 are shown in figure ITU1 
Asymptotically, as we would go higher on the m + -branch and rescale to restore unit 
mass configurations, the density would be completely pulled within the r = 1 region. 
The same was obtained by jTUj. 

7. Conclusions 

The Majumdar-Papapetrou formalism provides a good environment to study solutions 
to Einstein-Maxwell equations where matter is assumed to be described by electrically 
counterpoised dust (ECD). In this paper, we have obtained and analyzed regular and 
quasi black hole solutions stemming from the M-P formalism and obtained for diverse 
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Figure 9. The ratio <J g rr /gtt = \/A(r)/B(r) for the solutions shown in figure [8] 
(dashed lines) and for the ERN metric (dotted lines). 
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Figure 10. The density \/Ap(r) for the solutions in figure |H| The densities that 
effectively vanish at r ~ 1.1 correspond to m + — 100/rescaled solutions, while those 
that effectively vanish at r ~ 1.5 correspond to m + = 10/rescaled solutions. 



spherically symmetric ECD distributions by numerical integration of the nonlinear field 
equations. As an immediate consequence of nonlinearity, bifurcating solutions have been 
identified with respect to the amount of ADM mass allocated in the mass source term. 
Also an upper bound to the source strength has been found above which no solution 
exists. Although ECD distributions have assumed analytically different (spherically 
symmetric) forms, we have been able to reformulate the sources to treat them on equal 
footing. From this treatment we have been able to obtain regular solutions that come 
arbitrarily close to black hole solutions, the so called quasi black holes. Investigation 
of such bodies could be useful in investigations of interiors of black holes that still 
hide many unanswered questions. Bifurcation is not an unusual feature in gauge field 
theory fiQ] or in gravity fM\ 123] . It is encouraging that here we are able to show 
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that the bifurcation described in this work is not an artifact of a particular choice of 
charge/matter /energy density. Stability and bifurcation are closely related problems, 
so investigation of the above solutions with regard to stability is a natural extension of 
this work, bearing in mind the citation from Ref. ^3] '. . . stability analysis may be more 
expensive than the calculation of the solutions themselves.' 
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